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Renormalization of non-magnetic impurity potential by strong electron correlation is investigated 
in detail. We adopt the t-t'-t"-J model and consider mainly a 5-function impurity potential. The 
variational Monte Carlo method shows that impurity potential scattering matrix elements between 
Gutzwiller-projected quasi-particle excited states are as strongly renormalized as the hopping terms. 
Such renormalization is also seen by the Bogoliubov-de Gennes equation with an impurity, where 
the strong correlation is treated by a Gutzwiller mean-field theory with local renormalization factors 
and local chemical potentials. Namely, the 5-function potential is effectively weakened and broad- 
ened. We emphasize the importance of including the local chemical potential, which is paid little 
attention to in the literature, by physical consideration of the doping dependence of a local hole 
density. We also investigate effect of smooth impurity potential variation; the strong correlation 
yields anticorrelation between the gap energy and the coherence peak height simultaneously with 
large gap distribution, which is consistent with the experiments. 



I. INTRODUCTION 

Anderson's theorem tells us that the s-wave supercon- 
ductivity is insensitive to small potential scattering^. On 
the other hand, the d-wave superconductivity has zero su- 
perconducting gap in the nodal direction, and thus may 
be sensitive to disorder. However, experimental obser- 
vation of the high-temperature superconductivity, which 
many people are nowadays convinced has d-wave symme- 
try, seems robust against disorder^ 3 -^. For example, the 
high-temperature superconductors seem to exhibit more 
conventional behavior at higher hole doping rates where 
the systems are supposed to be more disordered. Further- 
more, the local density of states measured by the STM 6,7 
show clear V-shape at low energy that indicates the d- 
wave nodes are not much influenced by disorder. Theo- 
retically, it is proposed that this protection of V-shape is 
due to strong Coulomb repulsion between electron s 8 ' 9 ' 10 . 
Hence, detailed studies of effects of strong correlation for 
impurity scattering are necessary. 

In correlated systems, the model parameters are ef- 
fectively renormalized. For example, a hopping term is 
renormalized by a factor smaller than unity because hop- 
ping is more difficult in the presence of the double oc- 
cupancy prohibition. On the other hand, an exchange 
term is renormalized by a factor larger than unity be- 
cause each site is more often singly occupied. Then, the 
question is: how are impurity terms renormalized? In our 
previous paper—, we have presented preliminary results 
of the impurity renormalization. That is, local modula- 
tions of the hopping and the exchange term tend to be 
enhanced by the local renormalization factors, and impu- 
rity potential tends to be screened by the effective local 
chemical potentials originated from minimization of the 
total energy. In this paper, we investigate in more detail 
how impurity potential is renormalized, using the varia- 
tional Monte Carlo (VMC) method and the Gutzwiller- 
projected Bogoliubov-de Gennes (BdG) equatio n 11 ' 12 ' 13 , 
then discuss agreement and disagreement with the exper- 
iments. 



An extensive amount of literature has been devoted 
for the impact of impurities on the normal and super- 
conducting state of the cuprates. A detailed review dedi- 
cated to this subject was recently given by Alloul et alr^ 
We do not repeat the whole review here, but the theoreti- 
cal side may be summarized as follows: Suppose electrons 
in host metal interact with each other by the onsite re- 
pulsive Hubbard U terms, and let us put a non-magnetic 
impurity in it. Then, very strong impurity potential (uni- 
tary scatterer), such as a cavity in otherwise uniform sys- 
tems, tends to induce local magnetic moments near the 
impurity if U is large enough, whereas weak potential 
(Born scatterer) does not. 

However, the local moment formation in the super- 
conducting state seems relatively controversial. The mo- 
ments appear according to the theory based on the mean- 
field decoupling of the U term, e.g., by Chen and Ting 15 , 
and Harter et alr^. Although it can be a good approx- 
imation for small U, yet antiferromagnetic correlation is 
probably underestimated especially at large U because 
the superexchange process is not taken into account ex- 
plicitly. Considering that the local moments typically 
appear at sufficiently large U , it is critical to know in 
what range of U the theory is valid. An interesting con- 
trast is in the theory by Tsuchiura et al}— based on the 
one-site removed t-t'-J model, i.e., an effective U — > oo 
model. They adopted two different Gutzwiller approxi- 
mations that lead to two different BdG equations. One 
of them takes away the double occupancy prohibition in 
return for the uniform renormalization of model parame- 
ters. This calculation results in the local moment forma- 
tion. It also indicates that J may cause the appearance 
of the local moment even if U — because this effec- 
tive model is a "U = but finite J" system. In the 
other BdG equation of them, each local model parame- 
ter is dressed with an extended Gutzwiller renormaliza- 
tion factor that depends on the position. This calcula- 
tion in contrast results in the absence of the local mo- 
ments because electrons tend to avoid the impurity and 
the antiferromagnetism locally collapses. However, these 
local renormalization factors are, without local deriva- 
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tion, speculated from the previously derived formula for 
the uniform systemic, and may need to be verified in 
the future studies. Tsuchiura's work was criticized by 
Wang and Leei^, who applied an inhomogeneous slave- 
boson mean-field theory to essentially the same model. 
The result shows the local magnetic moments in the un- 
derdoped region. In the slave-boson mean-field theory, 
the double occupancy prohibition is relaxed by the sad- 
dle point approximation and only its average is satisfied, 
which probably reduces the influence of the spin-spin in- 
teraction effectively. To compensate it, Wang and Lee 
added a phenomcnological residual spin-spin interaction 
term, and the result depends on its magnitude. In ad- 
dition, Gabay et aL— recently obtained similar results. 
Liang and Lee^i applied the VMC and also concluded 
that local moments appear in the underdoped regime. 



The strong scatterers introduced above are modeled on 
in-plane impurities which substitute for Cu in the Cu02 
planes. Besides such unitary scatterers, all cuprate ma- 
terials are doped by out-of-plane ions that mostly occupy 
random positions in the crystal lattice or interstitial po- 
sitions. Such intrinsic impurities may be weak, but can 
be poorly screened by electrons around them because the 
cuprates are quasi-two-dimensional metal. In addition, 
these weak potentials are not expected to induce local 
magnetism. Influence of these "Born scattering poten- 
tials" on the local density of states was studied by, e.g., 
Wang and Lee^, and Nunner et al~2- In fact, however, 
the effect of electron correlation on them seems hardly 
discussed in the literature, with the exception of Garg et 
al£, Fukushima et al 10 , and Andersen and Hirschfeld 24 . 
That is a subject we would like to address in this paper. 



After we define our model in Sec. HH the renormaliza- 
tion of impurity-potential matrix elements is shown by 
the VMC calculation in Sec. lIIIl Then, the BdG equation 
based on the Gutzwiller approximation with local renor- 
malization factor s 11 ! 12 ' 13 are solved in Sec. IIV1 where 
we emphasize the importance of including local chemical 
potentials through the comparison with the method by 
Garg et al£ Note that these local chemical potentials are 
introduced for minimizing the total energy and are differ- 
ent from those used for the inhomogeneous slave-boson 
mean-field theor y 19 ' 20 ' 25 that are the Lagrange multipli- 
ers to enforce the average no-double-occupancy condi- 
tion. Rcnormalization of smoothly varying impurity po- 
tential is also presented to show anticorrelation between 
the gap energy and the coherence peak height compatible 
with large gap distribution, which is consistent with the 
experiments^^. 



II. MODEL 

We use the t-t'-t"-J model with an impurity term, 
namely, 



H = H tt 't" + Hj + -Himp, 



(1) 
(2) 

(3) 



where c ia (c ia ) is the creation (annihilation) operator of 
site i and spin a, and hi = c\ a Ci a . As the hopping, 
we take Uj — t,t',t", for the nearest, second, and third 
neighbors, respectively, and otherwise zero. The summa- 
tion in the J term is taken over every nearest-neighbor 
pair. The Gutzwiller projection operator P prohibits 
electron double occupancy at every site. 

In this paper, we focus on the renormalization of a sin- 
gle non-magnetic 5-function impurity potential located 
at i = 0, 
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kcr Ck'cT, 



(4) 



except for Sec. HVDi where we briefly discuss the effect 
of smoothly varying impurity potential. Here, Nl is the 
number of sites. 



III. VARIATIONAL MONTE CARLO 
CALCULATION FOR MATRIX ELEMENT 
RENORMALIZATION 

Here we calculate matrix elements of the impurity po- 
tential with respect to the uniform Gutzwiller-projected 
quasi-particle states of d-wave superconductors using the 
VMC method. 

Let us start from a uniform system without impurities. 
We assume that the ground state is well approximated 
by a Gutzwiller-projected d-wave superconducting state, 
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where 



Uk 



Ek 




Vk 




^ + A%, Afc = A v (cos k x ~ cos k y ), 

-2i v (cos k x + cos k y ) — At' v cos k x cos k y 
— 2i"(cos 2k x + cos 2k y ) — /i v . 
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iV e is the total number of electrons. The variational 
parameters A v , t' v , t", /i v are optimized so as to mini- 
mize the total energy. We also assume that the excited 
states are well represented by the projected quasi-hole 
wave functional, 
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|0). 



(6) 



Then, we are able to calculate the matrix elements and 
spectral weights using the excited quasi-hole wave func- 
tion with the ground state parameters. 

By switching on the impurity potential, these excited 
states should be mixed by the matrix elements, 

V k , k > = (fc T |4 T <v T |fc' T) + i-k' T l4|Cfc'il - k t). (7) 

We carry out VMC calculation for Vk,k> , and show that 
its renormalization by the strong correlation is similar to 
that of the total spectral weight, 
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FIG. 1: (Color online) Comparison of V ktk > /V k B k , s and 
Z k /Z k cs as functions of the hole concentration x in the case 
of (t',t",J)/t = (-0.3,0.2,0.3) optimized at each x. Each 
symbol represents transfer between the fc-points of the same 
symbols in the inset. 
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which is known to be strongly renormalizedSS. It is wor- 
thy to be noted that the BCS theory has predicted the 
matrix elements of the impurity potential and the total 



spectral weights are V k k , = u k u k ' — v k v k 



and Zl cs = 

1, respectively. On the other hand, the Gutzwiller ap- 
proximation (GA) yields renormalization, 
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where g t is the Gutzwiller renormalization factor for the 
hopping term. However, according to the conventional 
GA, Vk,k' is not renormalized because it originally comes 
from a particle number operator; on the contrary, a GA 
generalized for inhomogeneous systems yields renormal- 
ization as will be discussed in the next section. 

We plot V k ,k' /V k B k, s and Z k /Z BCS in Fig.Q]as functions 
of the hole concentration. Suppose the impurity potential 
is not too large. Then, the matrix elements V k _ k > perturb 
the system only if the excitation energies of \ka) and 
\k'a) are close. Therefore, here we plot Vk,k' connecting 
two symmetric reciprocal lattice points indicated in the 
inset of Fig. [TJ each pair of symbols in the inset refers to 
k and k' , which corresponds to the same symbols of V k ,k' 
and Z k — Zy in the plot. The variational parameters are 
optimized for each hole concentration. Note that V k ,k' is 
renormalized as strongly as Zk, and its renormalization 
factor is quite close to gt- Furthermore, Fig. [2] compares 
V k .k' of different bare parameters in the Hamiltonian. It 
suggests that the renormalization is insensitive to param- 
eters. 

These results may be understood as follows. Let us 
look at the Fourier transform of the (5-function impu- 
rity potential [Eq. (0)]. The sum of its diagonal terms 
(k = k') just slightly shifts the chemical potential. What 
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FIG. 2: (Color online) Comparison of V Kk > /Vffi of different 
bare parameters. Each symbol represents transfer between 
the fc-points of the same symbols in the inset. 



about the off-diagonal terms (k ^ fc')? If k were a site 
index, they would be renormalized as gt = 2x/(l + x) ac- 
cording to the Gutzwiller approximation, which is smaller 
than unity and is going to zero as x — > because it is 
more difficult to hop in the presence of the double oc- 
cupancy prohibition. Even in k space, if electrons are 
densely packed in the lattice, it must be similarly diffi- 
cult to hop from fc to a different fc'. Thus, the impurity 
potential should be renormalized by a factor similar to 
g t as we expected. 
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IV. GUTZWILLER-PROJECTED 
BOGOLIUBOV-DE GENNES EQUATION 

A. Renormalization by effective local chemical 
potentials 

We solve a BdG equation derived using the Gutzwillcr 
approximation with local renormalization factor a 11 ' 12 ' 13 
for non-magnetic systems. By requiring minimization of 
the total energy, the BdG Hamiltonian naturally con- 
tains effective local chemical potentials originating from 
the derivative of the local renormalization factors with re- 
spect to local particle densities. In the following, we show 
that the impurity potential is renormalized by those local 
chemical potentials. 

Let us assume that a good variational ground state 
can be represented in the form of P'\tp), where l^) rep- 
resents a wave function obtained later by solving a BdG 
equation. The operator P' contains a fugacity opera- 
tor to control the particle number as well as the original 
Gutzwillcr projection P. In the following, we use nota- 
tion, 
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for an arbitrary operator O. The built-in fugacities allow 
us to require conservation of the local electron densities, 
namely, 

(hi) = {hi) = n t . (11) 
Then, the Gutzwiller approximation yields 
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, (13) 



where Xi = 1 — rii. 

We consider non-magnetic systems where Xij = 

( c It c J't)o = ( c lx c il)o and A y = ( c j1 c *t)o are real 
numbers, and Ay = Ajj. Then, the total energy 
(H — [i J2i ^i) can be calculated by the GA as 



E GA = -4 dljtijXij - E 4 [ 2 ( 3 4 - 1)4 
+2(34- + 1) A?- + mrij] - fj, E th + V n , (14) 



where the summation of the kinetic-energy term is 
taken over every pair. Using Xij = J2a( c L c i<? + 

4 (T Ci<r)/4 and Aij = (c^c^ + ckc^ + c jiCi] + Cii c n )/A, 
the extremum condition of Eqa leads to a BdG 
equatio n 11 ' 13 represented by the mean-field Hamiltonian 
#BdG = J2(ij)XijdEoA/dXij + AijdEaA/dAij + 



J^ihidEGA/drit, namely, 

-ffBdG = - E g\jUjCi<r C ja ~ E J ( 3 9ij ~ l )Xi 
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E J(39ij + l)A y Ay - + tofa + F o"o- (15) 



Note that, in contrast to the conventional BdG equation, 
it contains the effective local chemical potential 
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where the summation of the J term is taken over the 
nearest neighbors of site i. By diagonalizing the BdG 
Hamiltonian, we obtain H B dG = J2 n =i ^n(7i n 7in + 
72„72n) + const., with E n > and 
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Then, |V) = n„7in72 Il |0), n % = 2£>? 
and A^ =E„ 
The inclusion of the local chemical potential /tj makes 
it harder to optimize the local mean-fields, and simple it- 
eration does not converge very well. Strategies to look for 
the minimum of the total energy E seem to work slightly 
better. We have solved the self-consistent equation for 
the systems of 24x24 sites with the periodic boundary 
condition. We set t' = -0.3t and t" = 0.2*. Then, J 
and /z are determined using the uniform system without 
the impurity so that x is the desired hole concentration 
as well as J(3g| J - + l)Ay = 0.3i. These values are fixed 
in solving the equation for the impurity systems, i.e., we 
neglect 0(1 /N£) shift of /z caused by the inclusion of the 
impurity. 

According to Eq. (jlip . the expectation value of hi is 
by definition not renormalized by any "g" factor as the 
hopping and the exchange term. However, as one can 
see in the BdG Hamiltonian in Eq. (fT5|) . the impurity 
potential can be compensated by fii- Therefore, we define 
a renormalized impurity potential by including difference 
of fj,i , namely, 



Vi = V Q S l0 - (Hi - /Zoo) . 



(18) 



Here, h<x> is /Zj— >oo and approximately equal to /z; °f t ne 
system without the impurity, which is nonzerc 29 . In the 
uniform system, however, one usually redefines ji+Hoo as 
the chemical potential, and ^oo does not explicitly appear 
in the calculation 2 ^. Figures [3] and 0] show the calculated 
renormalized impurity potential at the impurity site for 
various values of the bare potential Vq and the hole con- 
centration x > 0.05. We have also tried the systems with 
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FIG. 3: (Color online) The dots represent the renormalized 
impurity potential at the impurity site as a function of the 
bare impurity potential with the hole concentrations x — 0.05 
(blue), 0.125 (green), and 0.2 (red). The broken lines are gtVo- 
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FIG. 4: (Color online) The dots represent the renormalized 
impurity potential at the impurity site as a function of the 
hole concentration with the bare impurity potentials Vo = 
0.25t (red), It (green), and 1.5t (blue). The broken lines are 
9tV . 



x = 0.025, but were unable to reach the energy minimum 
possibly because there are a couple of meta-stable states. 

Note that Vo is strongly suppressed and is quite close 
to gtVo, where gt is the Gutzwiller renormalization fac- 
tor of the uniform system. These results agree with our 
VMC results in the previous section. Here, Vq deviates 
upward from gtVa when Vq is large near the half filling. 
This is possibly related to the position dependence of 
Vi. Originally, the impurity potential is nonzero only 
at the impurity site. However, after solving the BdG 
equation, the renormalized impurity potential distributes 
more broadly as shown in Fig. [5l Namely, as Vq becomes 
larger and x becomes smaller, the renormalization effect 
reduces at the impurity site whereas it broadens toward 
the neighboring sites. 

This broadening can be understood as follows: Basi- 
cally, energy loss by the impurity potential reduces elec- 
tron occupation at the impurity site. However, the hole 



prefers to move around to gain the kinetic energy. There- 
fore, to minimize the total energy, the <5-function impu- 
rity potential is broadened by /Xj. More explicitly, the 
contribution to fii from the kinetic energy contains a fac- 
tor 
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which behaves as ^/xjjxi near the half filling. Suppose 



the local hole densities behave as 



with 



some exponents 



6i, Sj. Then, dg\JdXi 



as 



x — ► 0. If Si ^ 5j, then /i, or fj,j diverges and the self- 
consistent condition is not satisfied. Therefore, <5; = 8 



which tends to make the hole distribution more uniform. 

Such extension of the impurity potential is also re- 
ported in the context of the unitary impurity potential 
by Poilblanc et al~ and in the weak-coupling context 
by Ziegler et alQ In addition, Bulut^ and Ohashi^ re- 
ported that the agreement between experimental data 
and the random-phase approximation is improved by 
adding phcnomenological extended range potential. 



0.3 
0.2 
0.1 



x=0.2. 

i 

\ 




; 0.125 ^ \ 

\\ 

0.05 \ V, 

1 


*^t= — 



(0,0) (1,0) 



(2,0) 

i 



(3,0) (4,0) (5,0) 



FIG. 5: (Color online) The position dependence of the renor- 
malized impurity potential for Vo = It. 

At half filling, the non-magnetic impurity potential 
should not affect the ground state because each site has 
to be occupied by one electron in any case. However, 
they do affect the ground-state energy. In the words of 
the BdG equation, the impurity potential is renormalized 
by fii in -ffBdG ; and not in Eqa- This is different from 
the well-known renormalization of the hopping and the 
exchange term described by gjj and gfy, these renormal- 
ization factors influence both the ground state and the 
ground-state energy. 



B. Local density of states 

The projected quasi-particle states P'jl n \ip) are ap- 
proximately orthogonal to each other—. We regard them 
as excited states and calculate the density of states 
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(DOS). Then, the local DOS (LDOS) is represented by 
N(r,u) = g>„ [Kl 2 % - E n ) + K\ 2 5(lo + E n )} . 

n 

(20) 

To obtain dense spectra, we use a supercell composed 
of 24x24 sites whose origin has the impurity, and this 
supercell is repeated so as to construct a superlattice of 
10 x 10 supercells with the periodic boundary condition 34 . 
Then, the Hamiltonian can be block-diagonalized by the 
Fourier transform with respect to the supercell indices, 
and the calculation of expectation values is reduced to an 
average over many quasi-twisted boundary conditions of 
the 24x24 site system. This supercell method is useful 
to obtain dense energy levels although it seems to over- 
estimate correlation functions between very distant sites 
if the supercell size is too small and the number of su- 
percells is too large. Except for this supercell boundary 
condition, the other conditions of the calculation are the 
same as those in Sec. lIV A"l Since spectra in finite systems 
are discrete, we replace each (5-function by the Gaussian 
distribution with the standard deviation SE = 0.02£ to 
obtain continuous DOS. 

Figure[Sfa) shows the calculated LDOS at the impurity 
site as well as its nearest, second, and third neighbors. 
The DOS of the uniform system (Vq — 0) is also plotted 
by dotted lines as a reference. Here, we have chosen Vq — 
It as the impurity potential, which is of the same order 
as the renormalized band width (8gtt). Nevertheless, it 
is well screened by /Lt« and the LDOS is not very site- 
dependent in agreement with results in Sec. IIV Al At 
the impurity site, the hole density is larger, and thus 
g* r is large. As a result the LDOS is also larger than 
those of the other sites. The peak at E ~ — 0.37t is the 
van Hove peak; the band renormalization by gt makes it 
closer to the Fermi surface. In fact, another van Hove 
peak appears in the positive bias at E ~ 0.37t. It is 
much weaker than the negative-bias van Hove peak, and 
only appears as a small shoulder of the the coherence 
peak in the resolution of Fig. [6^a). Some small portion 
of the original electron band around the van Hove peak 
is unoccupied in the mean-field superconducting ground 
state due to the electron-hole mixture. Even though that 
portion is small, the singular DOS enhances it to yield 
the positive-bias van Hove peak. 

For comparison, we also show in Fig. GDJb) the results 
of the system without the strong correlation, i.e, let us 
set g\- = <?|. = 1 and /i^ = in Eq. (fTS)) and solve the 
BdG equation. J and fi are determined in the same way 
as in the strongly correlated case, i.e., 4JAy = 0.3i for 
the uniform system. Here, we use supercells of 36 x 36 
sites to obtain dense spectra; if the same system size is 
used, energy level spacings are larger than those in the 
correlated systems due to the wider band width. The 
impurity potential is also Vb = It. Note that it is this 
time much smaller than the band width 8t. Nevertheless, 
comparing Figs.[SJa) and^b), it is clear that the system 
without strong correlation is much more disordered than 



the one with the correlation. 
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FIG. 6: (Color online) The LDOS for one impurity system 
with V — It and x — 1/8, at the impurity site and its nearest, 
second, and third neighbor, calculated by three different BdG 
equations: (a) with g\j, gL, g,i, (b) without correlation (g\j = 
gfj — 1, fj,i =0), (c) without the local chemical potential 
(with gjj, gfj, but \ii = 0). The dotted lines are the DOS in 
the uniform system (Vb = 0) as a reference. 

From Fig.[6jb), it is clear that the (5-function potential 
makes the LDOS asymmetrici 22 i 23 That is, some weight 
of the LDOS tends to move to the right at the impurity 
site, and to the left at the nearest neighbor site. The 
shift is large at high energy, and small at low energy. 
This asymmetry appears because the (5-function potential 
lifts the degeneracy of two quasi-particles: (i) a linear 
combination between an electron state at site 1 and a 
hole state at site 2, or (ii) a linear combination between 
an electron state at 2 and a hole state at 1. We discuss 
it more in detail using a two-site system in Appendix 
lAl With the strong correlation in Fig.^a), however, the 
asymmetry is less pronounced and the LDOS seems more 
symmetric. Most likely, the broadening of the impurity 
potential by fii results in weakening of the asymmetry. 

C. Comparison with other strongly-correlated BdG 
equations 

1. Importance of including the local chemical potential fii 

Garg et al. used a similar BdG equation^. Although 
the definition of the local renormalization factors is the 
same, they do not take into account the effective local 
chemical potential \ii that minimize the total energy. 
They have reported that the spatially averaged LDOS 
shows protected V-shapes at low energy, and stated that 
the correlated systems are less disordered. For the test- 
ing purpose, we have also solved their BdG equation and 
show the resultant LDOS in Fig. HJc). The parameters 
are the same as Fig. [6](a), and the uniform limits of these 
two BdG equations are identical. We have found that 
the LDOS before the spatial average is actually quite dis- 
ordered, especially at the impurity site, and that the V- 
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shaped LDOS appears only after the average. Therefore, 
it seems difficult to conclude that the LDOS in Fig. [S{c) 
is less disordered than the LDOS in Fig. [6] (b), in our 
opinion. 

A more serious problem of their method appears in 
the local hole density as a function of the bulk hole den- 
sity shown in Fig. [7] As already mentioned above, the 
non-magnetic impurity potential at half filling should not 
affect the ground state because each site has to be occu- 
pied by one electron in any case. Hence one can naturally 
imagine that the local hole density Xq at the impurity 
site approaches the bulk hole density X cLS X -» 0. How- 
ever, if ^ is not taken into account, x — x increases as 
x decreases. Therefore, it is questionable if the method 
by Garg et al. correctly captured the properties of the 
strongly correlated systems. 



0.15 



X 




0.05 



FIG. 7: The relative hole density xo — x at the impurity site 
as a function of the bulk hole density x solved for Vo = 0.25t 
by the two different BdG equations; without jj,^ and with fii. 



Capello et aL— solved the BdG equation with extended 
Gutzwiller renormalization factors, but without the lo- 
cal chemical potential. Although we have not duplicated 
their results, we speculate that they have similar prob- 
lems as Garg et al. because of the lack of the local chem- 
ical potential to minimize the total energy. 



2. Large U instead of the Gutzwiller projection 

Andersen and Hirschfcld 24 and some references therein 
used the BdG equation without the Gutzwiller renormal- 
ization factors (gt = g s = 1). They took into account 
the electron correlation by the Hubbard U term of the 
mean-field level. In that case, the screening effect similar 
to the one described in this paper can be obtained be- 
cause —Urii S plays a role of our local chemical potential 
\x%. Since the impurity potential reduces electron occupa- 
tion, the Coulomb energy loss is smaller there. Namely, 
summing up the potential and Coulomb energy loss, the 
total loss at the impurity site become less prominent by 
U. Note that, however, the strength of the screening 



depends on how one chooses U. A problem is that the 
mean-field decoupling of the U term may underestimate 
the exchange energy and overestimate the kinetic energy 
especially near the half filling; the ground state would be 
strongly influenced by them. 



D. Smoothly varying impurity potential 

As has been shown above, short-range impurity poten- 
tial is screened by fii ■ Accordingly, what remains must be 
the smooth variation in the potential. Here, we demon- 
strate the influence of the strong electron correlation on 
it using the Coulomb potential. 

Instead of Hi mp , we consider the Coulomb potential 



from randomly located off-plane impurities, namely, 



V c 



-J v /(r-r £ )2+d 2 ' 



(21) 
(22) 



where A^ mp is the number of impurities, rg is the posi- 
tion of ^-th impurity projected on the xy plane, and d is 
the off-plane distance. Here, we take Vc = U and d = 1 
and adjust /i to satisfy x = 1/8 simultaneously with the 
self-consistency condition. In addition, one supercell has 
20x20 sites with A^ mp = 12 impurities, and the same 
impurity configuration is repeated so as to construct a 
superlattice of 10x10 supercells. For simplicity, in deter- 
mining the Coulomb potential in Eq. (I22p . we use only 
one of the cells, i.e., the system of 20x20 sites with the 
periodic boundary condition. 

Figure M shows the LDOS at 4 sites A, B, C, and D, 
each of which has a different hole concentration. The 
LDOS in hole-rich regions (e.g., A) has high coherence 
peaks with low gap energy. In contrast, hole-poor re- 
gions (e.g., D) has low coherence peaks with high gap 
energy. Figure [5] compares the spatial dependence of the 
bare impurity potential V r — N^ 1 ^ r V r and the renor- 
malized impurity potential V r — N^ 1 ^2 r V r . Here, we 
have subtracted the spatial average of the potential. It is 
clear that the renormalized impurity potential has much 
smaller and smoother variation. 

Nunner et al£2- pointed out for the conventional BdG 
equation with smoothly varying potential that the LDOS 
at site r is similar to the DOS in the uniform system 
with shifted chemical potential jl = ji — V r . This argu- 
ment can be also applied to our system. However, the 
Gutzwiller renormalization factors make a difference in 
the relation among the hole density x r , the gap energy, 
and the height of the coherence peaks. In the systems 
without the Gutzwiller projection, is determined by 
the DOS near the Fermi level. Namely, Ay takes max- 
imum when the chemical potential is at the van Hove 
singularity. Accordingly, when A^ is large, the gap en- 
ergy and the peak height are also large. In contrast, 
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FIG. 8: (Color online) (Left) The LDOS at 4 sites A, B, 
C and D indicated in the right figure, for the system with 
smoothly varying potential with x = 1/8, Vc = U, d = 1, 
]Vi m p = 12, and supercell size 20x20 sites, solved by the 
Gutzwiller-projected BdG equation. (Right) The hole den- 
sity distribution in a supercell. The white dots are positions 
of the impurities. 
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FIG. 9: (Color online) Spatial dependence of the bare impu- 
rity potential Vj — N£ £\ Vj and the renormalized impurity 
potential Vi — £\ K along a diagonal line (from the left 
bottom to the right top in the hole density plot of Fig. [8} ■ 



with the Gutzwiller renormalization factors, as i — » 0, 
(i) the band shrinks by and thus the DOS near the 
Fermi level increases, and (ii) the pairing interaction is 
enhanced by gr?.. Because of (i) and (ii), Ay and the gap 
energy monotonically increases as x — > 0. Furthermore, 
the LDOS contains an extra factor g* i: and thus the peak 
height decreases as x decreases. 

Note that this anticorrelation between the gap energy 
and the peak height as well as the large gap distribution is 
consistent with the experiment er' 26 . However, according 
to our calculation, the large spatial variation in the LDOS 
accompanies large variation in the hole density, which has 
yet to be verified by the experiments. 



V. SUMMARY AND DISCUSSION 

In this paper, we have investigated the renormaliza- 
tion of the impurity potential by the strong correlation. 
The VMC calculation has shown that impurity potential 
scattering matrix elements between Gutzwiller-projected 
quasi-particle excited states are as strongly renormalized 
as the hopping terms. It may be understood by the 
Fourier transform of the 5-function impurity potential 
having a form of the hopping term in k space. The (real- 
space) hopping term is known to be renormalized by a 
factor less than unity because it is more difficult to hop in 
the presence of the double occupancy prohibition. Even 
in k space, if electrons are densely packed in the lattice, 
it must be similarly difficult to hop from k to kl ^ k, then 
the impurity potential and the hopping term should be 
renormalized similarly. 

Such reduction in the impurity potential is also seen by 
the BdG equation with local Gutzwiller renormalization 
factors. Near the half filling of the strongly correlated 
systems, the influence of the non-magnetic impurity po- 
tential on the ground state is small because each site has 
to be occupied by almost one electron in any case. How- 
ever, the impurity potential does affect the ground-state 
energy. Such properties appear by taking into account 
effective local chemical potential, which is paid little at- 
tention to in the literature. In addition, the local chem- 
ical potential effectively broadens the impurity potential 
because holes prefer to move around to gain the kinetic 
energy. Effect of smoothly varying impurity potential 
has been briefly discussed. It shows large gap distribu- 
tion. If the Gutzwiller renormalization factors are taken 
into account, the gap energy and the peak height are 
anticorrelated. These properties are consistent with the 
experiments 6 -^ 6 -. 

In fair comparison, there are also some disagreements 
with the experiments. According to our results, short- 
range non-magnetic potential variations are reduced, 
thus the system is more uniform, and accordingly the d- 
wave superconductivity can be robust. However, in the 
experiments, the system seems quite disordered but still 
the d-wave is robust. Such short-range disorder may be 
introduced by spatial modulation of t%4 or J, which can 
be enhanced by the strong correlation^ -, or by magnetic 
impurities, or by the electron-lattice interaction^ 6 .. 

In addition, although the smooth impurity potential 
variation yields anticorrelation between the gap energy 
and the peak height, it does not show the almost spa- 
tially independent V-shaped LDOS at low energy seen 
often in the experiments, which can be explained instead 
in the case of the rapid potential variational. In our 
previous paper—, we have discussed the LDOS of stripe 
states, where Ay contains two components; one is spa- 
tially uniform and the other is oscillating, typically with 
wave number q = 7r/4 or n/2. Then, the V-shaped gap is 
determined by the uniform component, and the oscillat- 
ing component influences it little. As a result, the linear 
slope of V-shape is robust (it does not have much spa- 
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tial dependence) . This is rather counterintuitive because 
one may think that the local gap could be determined 
by Ay. However, the superconducting gap is not de- 
termined by such local properties, but by "coherence" , 
i.e., spatial dependence of Ay. Let us recall that, in 
the case of the zero-momentum pairing as the conven- 
tional BCS theory, the spin-up electron band couples 
with the spin-down hole band (upside-down spin-down 
electron band); these bands intersect at the Fermi level, 
and a gap opens by a nonzero superconducting order pa- 
rameter. On the other hand, the oscillating components 
contain only pairing with nonzero center-of-mass momen- 
tum. Then, the spin-up electron band couples with the 
g-shifted (or multiples-of-g shifted) spin-down hole band. 
The point is, these bands typically intersect not at the 
Fermi level except for limited points. In such cases, the 
V-shape of the LDOS is not affected a lot. 

Similarly, oscillations of the variational parameters 
other than Ay with wave number q mix "the bare band" 
and "the bands shifted by multiples of q" . Such terms 
change the band structure especially near band intersec- 
tions. However, this change in the band structure is not 
related to the superconductivity, at least in the mean- 
field approximation. The superconducting DOS depends 
on where one puts the chemical potential, but general 
properties near the V-shaped LDOS is determined by 
Ay- 

The smooth impurity potential variation considered in 
Sec. lIVDl is similar to the stripe states with q ~ 0. Since 
q is small, the oscillating components of Ay has effect 
similar to the uniform component q = 0. Namely, it 
affects every state at the Fermi level. However, the gap 
has now spatial dependence because q is not completely 
zero. In this case, the local gap is determined by Ay-; 
the LDOS at sites i is similar to the DOS in the uniform 
system with A x = -A v = A^ i+£ , x = Xij> and A = 
(i — Vi. Therefore, there is nothing like the "uniform 
component" in the argument of the stripe state, and thus 
the linear slope of V-shape in the LDOS is not robust 
anymore. This issue will be addressed in the future. 



Acknowledgments 



The authors thank C.-M. Ho for stimulating discus- 
sions. This work was supported by the National Sci- 
ence Council in Taiwan with Grant no.95-2112-M-001- 
061-MY3. Part of the numerical calculation was per- 
formed in the IBM Cluster 1350, the Formosa II Cluster 
and the Triton Cluster in the National Center for High- 
performance Computing in Taiwan, and the PC Farm II 
of Academia Sinica Computing Center in Taiwan. 



APPENDIX A: LDOS ASYMMETRY CAUSED 
BY (5-FUNCTION POTENTIAL 

As shown in Fig. [6l the (5-function potential causes 
asymmetry in the LDOS near the impurity site. Such 
asymmetry may be important because some STM exper- 
iments analyze data by taking the ratio between intensi- 
ties of positive and negative bias^: the symmetric part of 
the LDOS is canceled and its asymmetric part remains. 
To explain the origin of the asymmetry, let us diagonalizc 
the BdG Hamiltonian of a simple two-site problem ana- 
lytically in the following. We know that it is not realistic 
to apply the mean-field approximation to two-site prob- 
lems. However, it provides us some insight into numerical 
solutions in the bulk systems. 

Suppose the potentials at sites 1 and 2 are V (> 0) 
and —V, respectively. By putting the states in the order 
of electrons 1 f , 2 f , and holes 1 j, 2 j, the BdG Hamil- 
tonian matrix is written as 



Ho 



V -t A\ 
—t -V A 
A -V t 
AO t Vj 



(Al) 



Let us start from a simple case of t = and assume 

V < A. Then, A mixes 1 f electron and 2 J, hole with 
the equal weights. We call the linear combination of 
them with positive and negative energy the quasi-particle 
A + and A~, respectively. Similarly, B ± denotes the lin- 
ear combination of 2 f electron and 1 j hole with posi- 
tive/negative energy. 

If V = 0, A and B are degenerate. Finite but small 

V causes energy difference between A and B, namely, 
E A ± = ±A + V, E B ± = ±A - V. In the ground state, 
A - and B~ are occupied, but A + and B + are unoccu- 
pied. First we focus on If . Then, the electron addition 
spectrum is at E = E A + , the removal is at E = E A - . 
Namely, the finite V just shifts the whole spectra by V, 
and the local spectra are not symmetric around E = 0. 
On the other hand, the spectra of site 2 shift by —V. As 
a result, the spectra summed over sites 1 and 2 are still 
symmetric. As for If, the addition is at —Er- = E A +, 
and the removal is at — E B + = E A - (these negative signs 
originate from the treatment of If as a hole), which are 
identical to the spectra of If as expected. 

Next, let us consider f / 0. Then, t mixes A + 
and B _ , and the eigenstates are linear combinations 
of them whose eigenenergies are ±[i 2 + (A + V) 2 ] 1 / 2 . 
The eigenenergies of superposition of A - and B + are 
±[t 2 + (A - V) 2 ] 1 / 2 . As t increases, the LDOS asymme- 
try becomes weaker, but it never disappears. 

Back to the bulk systems, maybe we can physically 
explain the asymmetry as follows. Since J term causes 
pairing, a Cooper pair is formed more or less between 
nearest neighbors (site 1,2) when a snapshot is taken. 
This Cooper pair is a resonance of "states where sites 1 
and 2 are occupied by electrons" and "states where both 
are occupied by holes" . Destruction of the pair leaves 
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a quasi-particle. There are two possibilities for it; (i) 
an "electron at site 1 and hole at site 2" , and (ii) an 
"electron at 2 and hole at 1" . The ^-function potential 
lifts the degeneracy of these quasi-particles. That re- 
sults in the asymmetry in the local spectra although the 
bulk spectra are symmetric as in the two-site problem 
above. In the bulk system, the spectra are continuous, 



and the spectral shift by the impurity potential is larger 
at high energy than at low energy because (i) the near- 
nodal quasi-particles at low energy are less influenced by 
the impurity potential because there are not many states 
to mix with, and (ii) the shift is caused by A and it is 
smaller at lower energy. 
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